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Abstract 

The 2D Hubbard model is extended by placing 25* orbitals at each lattice site 
and studied in a systematic 1/S expansion. The 1/S results for the magnetic 
susceptibilitity and the spectra of spin-wave excitations at half filling are con- 
sistent with the large S calculations for the Heisenberg antiferromagnet. The 
1/5 corrections to the fermionic spectrum lift the degeneracy along the edge 
of the magnetic Brillouin zone yielding minima at (ib7r/2, ib7r/2). Relation to 
previous papers on the subject is discussed. 



1 



Typeset using REVT^ 



I. INTRODUCTION 



The 2D Hubbard model and the closely related t — J model have been extensively used 
recently in the context of high temperature superconductivity as the simplest models which 
capture at least some of the exciting physics of high— Tc superconductors. At half-filling, 
the large-f/ Hubbard model reduces to the nearest-neighbor Heisenberg antiferromagnet [|T| , 
which perfectly describes the properties of pure La2CuOi while at sufficiently strong 

doping, it describes a strongly correlated metal and possibly a superconductor [TyjTTj]. 

One of the crucial issues for lightly doped antiferromagnets is the role of quantum fluc- 
tuations. In La2-xSrrcCu04, they destroy long-range order in the ground state already at 
a very small doping f^. There is therefore a need for a systematic accounting of quantum 
fluctuations in the Hubbard model. For ordered Heisenberg antiferromagnets, a customary 
tool for this is the 1/5* expansion. This expansion gives rather fast convergence and first 
order 1/S results are quite accurate even for S = 1/2 0,0] • 

In this paper, we report the results of the large-S" studies of quantum fiuctuations in the 
Hubbard model at half-filling. We view the results as a basis for future considerations of 
quantum fiuctuations in doped antiferromagnets. 

To perform 1/5* studies, we need to extend the Hubbard model to higher spins. A 
straightforward way to do this is to place 25' equivalent orbitals at each lattice site |T^,|T3[. 
The hopping term remains the same as in the conventional Hubbard model, only the hopping 
integral t now scales with S. At the same time, the Un^rii term is replaced by the Hund 
rule interaction term 

Hint = ^ I] ^Ln^lpJ^iPn^jc^n (l) 
n,a,i 

Here a = 1, 2 is the spin index, and i is the orbital index which runs over « = 1, 2.... 25*. For 
S = 1/2, the summation over spin index is absent and (|l]) reduces to the conventional on-site 
repulsion term. For S > 1/2, the coupling term ensures (at large U) a parallel alignment of 
the spins of electrons at all orbitals, i.e., total spin S at each site, which forms a natural basis 
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for 1/S expansion, and also ensures an equal number of electrons at each site. Fluctuations 
of this number cost energy of the order of US ^ t. The interaction term is simply related 
to the total spin per site as Hint = —UJ2n{Sn + {J^n — 2S')^)/4 where Afn = i^lani'ian is 
the fermion density, and Sn = X^j Si„ = \ip\an(^aP'^ii3n is the total spin on a given site. The 
fluctuations of the square of the total density are irrelevant at large [/, and we have checked 
explicitly that one can equially well use either the eq. (|l]) or —UY^n for the interaction 
term. The former choice is however advantageous for computational purposes. 
In the momentum space, the extended Hubbard Hamiltonian takes the form 

H = '^eka]ak'^jak + X! '^Ia,k''^ll3,-k'+g'^il3-k+q0'ja,k (2) 

k k,k',q 

Here = —2t{coskj. + cosky), N is the total number of lattice sites, and the momentum 
summations extend over the first Brillouin zone — vr < k^, ky < n. 

The large-S" expansion for magnetization and spin susceptibility of the Hubbard model 



has been previously studied in |]I5| . Although a complete agreement with the similar expan- 
sion in the Heisenberg model was obtained, we believe that several relevant contributions 
were missed, and the agreement is thus a bit fortuitous. Below we will point out the discrep- 
ancies with Ref. |T^ explicitly. The semiclassical expansion results for the magnetization 
were also reproduced in ||13| by means of a expansion, where z is the number of nearest 
neighbors. 



II. MEAN-FIELD THEORY 

We now briefly outline the results of mean-field {S = 00) studies of the Hubbard 
model P,p|, p!2| , p!3| , p!6| . At the mean- field level, one introduces the condensate of Sz{q = Q) 
where Q = (tt, tt), and uses it to decouple the interaction term in (|^). The resulting quadratic 
Hamiltonian then takes the form 

-f^MF = ^kAak(^i<yk ^"^ Aa,k+Q^al3^iPk, (3) 

k ^ A; 
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where the summations over momentum here and below are hmited to the magnetic zone, 
i.e. half the first Brillouin zone. Note that < 5*^ > is the exact value of the condensate 
which includes all zero-point fluctuations. This quadratic Hamiltonian can be diagonalized 
by means of a Bogolubov transformation 

Ciak = Ukttiak + VkO'ai30'il3,k+Q 

diak = VkO-iak — UkO'af30"i(3,k+Q (4) 

where c represents conduction, and d — valence electrons. The transformation explicitly 
reads 

-G(-t))" -G(-t))"^ 

where = (e| + A^)^/^ and A = . The diagonalized Hamiltonian then takes the 

form 

H = Y, Ek{c\^^Ci^k - dl^kdiak)- (6) 

k 

The value of the gap A is obtained from the self-consistency condition of the diagonalization 
procedure 

- = y—. (7) 

U k^k ^ ^ 

The bosonic spin-wave excitations appear in the theory as poles of the total transverse 
susceptibility 

q'; u) = i jdt < TS^{t)S-_,m > e^"*. (8) 

To zeroth order in l/S, the total static susceptibility is given by a sum of bubble diagrams 
in Fig. 0: 



Xo 



) jY X! ^ EkEk+q ) -Efc + Ek+q ^ 



This simple form familiar from paramagnon theories exists, however, only for the static 



susceptibility. The total dynamical susceptibility is the solution of a 2 x 2 problem |T6 
because the antiferromagnetic ordering doubles the unit cell, and Xo i^^''^^) is nonzero 
either when q = q' or when q = q' + Q {in the latter case Xo'iQ^I + Q,uj) ~ u). Also 
notice the overall factor of S in Xo~ ■ It comes from the summation over the orbital index 
in each bubble. Because of this factor, Uxo^ is 0(1), and this makes the total transverse 
susceptibility very different from the bare one. 

For large A/t, one can expand in ( |10|) and simphfy (|]) to 

x^'iq) = — 7^ — r (11) 

^ 2J(l + 7,) ^ ^ 

where 7^ = {coskrc + cosky)/2 and J = 4t'^ /{2S)'^U. As q approaches the ordering momentum 
Q, 7g —1, and static susceptibility diverges as it, indeed, should because of the Golstone 
theorem. 

III. 1/5 EXPANSION FOR THE STATIC SUSCEPTIBILITY 

We now proceed to the study of the fluctuation corrections to the static susceptibihty. 
These corrections are related to the residual interaction between fermions. In the transverse 
channel there are two terms — the direct interaction of two fermions with opposite spins, 
and interaction between fermions of opposite spins mediated by magnetic fluctuations, see 
Fig.||. The effective Hamiltonian for fermion-magnon coupling was derived by D. Frenkel 



and one of us [|T^]; here we quote only the result: 



q) 

k,q 

+4akdif3,k+A<} '^^dik, q) + 4afeCi/3,fc+ge]3g ^dc{.k, q) + H.c. ) (12) 

Here Caq is the boson annihilation operator with polarization a, and the vertex functions 
are given by 

^cc4d{k, q) = ±(efc + ek+q)r]g + (e^ - ek+q)f] /2S 



$cd,dc(/c, q) = U 
and r]q and r/^ are given by 



1- o.r h^^T 1- 



8A2 " \ 8A2 



(13) 



In ([T3| ) we neglected all terms of the order J{t/US)"^ ^ J with m > 2 which is consistent 
with our intention to calculate corrections on the scale of magnetic interaction J. Further, 
unlike fermionic energies, the vertex functions do not have the overall factor of S. The 
second-order self-energy corrections to the fermionic propagators then will contain one power 
of S in the numerator due to the summation over the orbitals, and one power of 5* in the 
denominator which comes from Ek oc S. Total corrections will therefore be 0(1) i.e. will 
have 1/S factor compared to the mean-field inverse propagators which scale as 0{S). 

To calculate the spin susceptibility we will need the renormalized self-consistency condi- 
tion, i.e. expression for the order parameter < Sz > itl terms ofA = U<Sz>. It should 
ensure the cancellation of the leading 0(1) term in the denominator of • At S = 00, this 
condition was given by eqn. (^. For 1/5* corrections, however, we need the full expression 
for < Sz >. Substituting (|^) into the formula for the z component of the spin, we obtain 
after a simple algebra 

< Sz>= ;§-(< 4akdiak > - < cl^kCiak >) + ^(< cl^k^iak > + < 4^^Ciak >) (15) 

At 5 = oo, the only nonzero average is < d\^j^diak >= 1, and using A = U < 5*2 >, we 
return to (|^). However, at finite S, one has to calculate pair averages with the self-energy 
corrections. These corrections include terms which transform a conduction fermion into a 
valence one and vice versa after emitting or absorbing a spin wave. Accordingly, all pair 
averages, including < c\^j^diak > will be different from zero at finite S. The density diagrams 
are presented in Fig ^. A straightforward calculation using (|12D and ([13|) yields 
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where Qq = AJS Jl — 7^ is the bare spin-wave frequency, and we defined 




(17) 



Retaining only terms of the order of U we arrive at a well-known result for the Heisenberg 
model 



We are now in a position to compute 1/S corrections to the static transverse suscep- 
tibility. Simple considerations show that to the first order in l/S*, the RPA approach is 
still exact, but one has to include all 1/S corrections within a bubble and also take into 
account the renormalization of the relation between A and U . A special care should be 
given in selecting the diagrams as the valence-conduction bubble itself has the order oil/U. 
However this leading contribution gets cancelled in the denominator of (y), and the mo- 
mentum dependence of is given by the subleading terms in Xo~ which have the order 
of J/U"^. The cancellation of the 0{U) terms should clearly survive in l/S* expansion, and 
we, therefore, have to keep all 1/S diagrams to the order J/U'^S. Further, one should also 
consider the effect of including the direct longitudinal fermion-fermion interaction. Simple 
power counting arguments show that this does not give rise to an extra 1/S smallness be- 
cause each time one includes an extra interaction, one also has to sum over the orbitals of 
intermediate fermions. We found however that in the large-f/ limit each inclusion of more 
than one triplet ineraction does produce a smallness, not in l/S* but in (t/US)'^, because the 
frequency integration selects only those combinations of conduction and valence fermions 
for which the actual longitudinal interaction becomes small after being dressed by the Bo- 
golubov coefficients. Altogether, we found 14 relevant diagrams. They are presented in 
Fig ^ - The diagrams in Fig.^ contain no spin wave propagator, but two direct fermion- 
fermion interactions. The diagrams in Fig.^ contain self-energy corrections related to the 
exchange of spin-waves and no direct fermion-fermion interactions. The diagrams in Fig.|^ 




(18) 
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contain one spin-wave propogator and one direct interaction between fermions with parallel 
spins. Finally, the diagrams in Fig.^ contain one spin-wave propagator and two interactions 



in the triplet channel. The latter diagrams were omitted in Ref together with several 
0{J/U'^S) terms in other diagrams. We found that altogether, the terms omitted in |15 



cancell each other, but we do not believe that this could be anticipated in advance. Note 
that inclusion of self-energy corrections with the direct transverse interaction would lead to 
double-counting, as these corrections have already been taken into account in the mean field 
diagonalization. 

Calculation of the diagrams is tedious but straightforward. We list the expressions in the 
Appendix and quote here only the final result. We found that the modified self-consistency 
condition ensures an exact cancellation of the denominator in (j^) at g = Q in accordance 
with the Goldstone theorem. The functional form of x'^ is also unchanged compared to the 
mean-field result, but quantum fluctuations give rise to an overall renormalization factor. 



Z^. We have 



y2 
1 

As expected, this result coincides with the known 1/ S expression for the static susceptibility 
of the Heisenberg antiferromagnet 

Finally, for completeness, we will also obtain the 1/5* expansion result for the spin- wave 
velocity, c. From hydrodynamical considerations fl^, we know that (? = Ps/x± ? where x± 
is the static susceptibility 

X± = lx^-iq = 0,u = 0) = ^, (20) 
and ps is the spin stiffness defined by 

x::iq^Q)= (21) 

Ps{q-Qy 

Using the result for Sz and expanding (|19|) near q = Q, we find 



o-JS'Z- Z -l-M + £-i_ J_vi^ii— ^ + — (22) 



Combining (pOD and p^), we obtain 

c^ = 2J^Zl- ^.^ = |^ = l + |(5-A) = l + A^[l_yrT^] (23) 

which again coincides with the known spin- wave result. 

IV. SELF-ENERGY CORRECTIONS TO THE FERMIONIC SPECTRUM 

We now turn to the calculation of self-energy corrections to the single-particle excitation 
energy, Ek- In the mean-field approximation we had Ek = y^A^ + e|. At large U, it 
can be simplified to -E^ = A + AJS{coskx + cos/c^)^. These forms for E^ imply that the 
minimum of the fermion energy coincides with the whole edge of the magnetic Brillouin zone, 
\kj. ± ky\ = 71. This degeneracy in the position of the minima of Ej^ was responsible for the 
singular behavior of the Hubbard model at infinitesimally small doping ^j. We will show, 
however, that this degeneracy does not survive perturbative 1/S corrections, and the actual 
minima of the quasiparticle energy are located only at the particular points (±7r/2, ±7r/2) 
in the Brillouin zone. This is consistent with the results of other approaches pT|-[25|]. 

Consider for definiteness the corrections to the energy of a valence fermion. As before, 
the dominant self-energy corrections at large U are those which include the exchange of 
transverse magnetic fluctuations. To first order in 1/S, there are two such diagrams (Fig.^. 
Performing the frequency integration in these diagrams, we obtain for the Green function 
of a valence fermion 

G-\k, co)=co + E,-§J: ( ^ ^--^ ^ ^ + -] (24) 

We now consider the two self-energy terms separately. The first term in (p^) corresponds to 
a process when valence fermion transforms into a conduction one after emitting a spin wave. 
The matrix element for such process and the energy denominator scale as U and U S, corre- 
spondingly. As a result, the first term yields the momentum-independent renormalizations 
of the gap, A A(l + A/S), and of the wavefunction residue, Z = 1 — A/2S. The latter is 
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in fact related to the renormalization of the sublattice magnetization |T3[ < Sz >= Z — 1/2. 
This first self-energy term also has momentum-dependent contributions due to the k depen- 
dence in Ek+q, but they scale as J and, as we will see, are comletely overshadowed by the 
momentum-dependent contributions from the second self-energy term. 

The second self-energy term in (|2^) describes a process which includes only valence 
fermions. The vertex function for such process scales as t/S so intuitively one may conclude 
that the second term is less relevant than the first one. However, at resonance, u = —Ek, 
and the leading 0(A) terms in the denominator get cancelled. The remaining terms scale as 
JS, so that the total self-energy correction behaves as 0{A/S). More important, however, 
is that the k— dependent self-energy terms are not anymore suppressed by a factor J/U 
and, therefore, may substantially modify the shape of the excitation spectrum. 

The general structure of the 1/S self-energy correction is rather involved, so below we 
will restrict calculations to the dispersion relation along the lines \kx ± ky\ = tt where the 
mean-field energy of valence fermions has a degenerate maximum. Simplifying eq.(^^ at 
7fc = —1 and US/t ^ 1 we obtain for the pole of the Green function 

Et''^A(,-Um. (25) 



S 2S 



where 



2 _ 1 + a/1 - 7, 



/^W = ]^E , , 1 H \ - (26) 



It is not difficult to see that f3{k) does contain some momentum dependence, and, therefore, 
the excitation energy acquires dispersion along the edge of the magnetic Brillouin zone. We 
evaluated the integral (pB]) numerically and found that j3{k) has four equivalent minima at 
(±7r/2, ±7r/2). Near each minimum, one can expand in kx,y = ±7r/2 — k^^y and obtain the 
dispersion relation typical for the anisotropic 2D Fermi gas: 

Here A differs from A due to momentum-independent self-energy corrections, and we in- 
troduced Pj_ = {kx ± ky)/\/2, p\\ = {kx =F ky)/^ (the upper sign is for (7r/2,7r/2) and 
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(— 7r/2, —1^/2) and the lower sign is for the other two minima). One of the effective masses 
is finite aheady at the mean-field level, and we easily obtain from E^^°°: 

mx = ^ (1 + 0(1/5)) (28) 

On the contrary, mw is infinite at the mean-field level, and acquires a finite value only due 
to 1/5* corrections. The numerical evaluation of this mass yields 

2S 



mil 



(29) 



0.084A 

For the ratio of the effctive masses we then obtain to the leading order in 1/S: m\\/m± ~ 
188(JS'/f/). Notice the surprisingly large numerical coefficient. 

Care has to be taken in applying the 1/S results for the masses to the physical case 
of S* = 1/2. We already mentioned that in resonance, the diagram with only valence 
fermions has small energy denominator 0{JS) because the leading A term in -Efe+g gets 
cancelled by the external frequency. In this situation, vertex and self-energy corrections to 
the internal fermion line are important as they contain large factors ofU/J. A simple analysis 
similar to that in shows that self-energy and vertex corrections scale as 0{U/JS) and 



0{U/JS)^ respectively. Neglecting these corrections is therefore legitimate only if U/JS ^ 
1 and, in a rigorous sense, our 1/S expansion for the fermionic energy is valid only if 



S ^ U / J . It has, however, been argued several times [|l^,|l^,^] that at large [/, 
the actual (renormalized) interaction between fermions and spin-waves (which includes the 
wavefunction renormalization factor of fermions) retains its functional form (at least near 
magnon momentum q = (7r,7r)), but has the same order of magnitude as the bandwidth. 
The argument here is that for f/g// ~ JS, self-energy and vertex corrections do not generate 
any new energy scales besides JS. In this situation, the actual expansion parameter for our 
studies is 0(1), and we can expect that the lowest-order approximation we are using still 
yields at least qualitatively correct results. 

We conclude this section by a brief comparison with other approaches. Vignale and 
Hedayati PT| considered the same perturbative corrections as we did, but resorted to self- 



consistent rather than perturbation analysis. They required that all corrections to the bare 
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disperson of quasiparticles be on the scale of JS, and solved self-consistently Dyson equation 
for the quasiparticle energy. They found the band minimum at (7r/2,7r/2), which agrees 
with our finding. We, however, did find the corrections of the order 0{A/S) in perturbative 



calculations. Boninsegni and Manousakis ||2^ used variational Monte-Carlo to compute the 
hole dispersion in t — J model. They also obtained band minima at (±7r/2, ±7r/2). The 
quantitative comparison with our results is, however, difficult because of the relatively small 
values oit/J < 0.5 considered in The minima at (±7r/2, ±7r/2) in the t — J model 
were also found in the perturbative calculations in the small t/J limit Trugman ||24 



and Sachdev ||25[ performed variational calculations for the Hubbard and t — J models 
respectively and also found the minima at (±7r/2, ±7r/2) at large t/J. We compared the 
bandwidth along k^ + ky = n and found good agreement between Sachdev's results and ours. 
Thus at t/J = 2, variational and our 1/S analysis yield AE = E{tt, 0) — E{7t/2, 7t/2) ^ 0.2t, 
and 0.24t correspondingly. 



V. CONCLUSIONS 

To summarize, in this paper we considered a large S extension of the Hubbard model and 
calculated in l/S" expansion the leading quantum corrections to the static spin susceptibility, 
spin-wave velocity and hole dispersion at half-filling. We found that the mean-field degen- 
eracy of the quasiparticle spectrum does not survive in l/S* calculations, and the actual 
hole dispersion has minima at (±7r/2, ±7r/2). Our results for magnetic parameters agree 
(as they should) with the 1/S calculations in the Heisenberg model. Indeed, at half-filling, 
the use of Hubbard model for magnetic calculations is not the easiest way to arrive at the 
final result. However, the method we presented here has a substantial advantage over other 
techniques in that it can be straightforwardly extended to doped antiferromagnets where 
quantum fiuctuations are much more relevant than at half filling. The work along these lines 
is now in progress. 
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VII. APPENDIX 

Here we give the expressions for the diagrams in Fig. § - ^. The two diagrams on 
Fig. ^ are obviously irrelevant to magnetism, and, accordingly, cancel each other out giving 
contributions of i^^^.Two of the rest of the diagrams, namely those of Fig. |^ and |^a, 
contribute to order 1/US, while the rest contribute to order J/U'^S. The explicit expressions 



for the first two diagrams are 



Xla 



-E- |i 

JV2 I 



4A2 



(30) 



for the diagram |a, and 
1 



X2a 



E 



[l + l,) A- ^ B 



(31) 



for the diagram Here all the sums are over a half of the Brilloiun zone, and constants A 
and B were definied in ([T7|) . Using the self-consistency condition (|16|), the total contribution 
from these diagrams can be immediately simplified to 



Xla + X2a 



1 

u 



E 

k 



2S 

Ek 



2J 



{2A - B) 



2J 



(1 + 7.) A 



(32) 



Care also should be given to the calculation of the bare bubble, because the relation between 



A and U contains 1/5* terms. From (]T0|) we have 

2S 2J , ^ fUS\^ ^ 2S 



Xo = Xq 



E 



Eh 



f/2 



Eh 



2J 
IP 



7,J 



3A\ 



(33) 



Assembling eqn (|32D and (0), we obtain 
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1/27 \ 27 47 

Xo + Xi + X2 = jj (l--^(l + 7.)) - jj^{2A-B) - ^(1 + 7,)A (34) 

The expression for other diagrams are as follows 

Xlb= 7^ ^ H V^I-tI; X1c + XW= (l + 7g) (^-5); Xle = (l + 7g) ^ 

A; 

2J 2J 2J 

X2fe = - JJ^ Iq - B)- X2c = Iq B] X2d = (1 + Iq) B] 

^2e = - ^ E /i^tI; X2f = - (1 + Iq) A. (35) 

Finally, the diagrams in Fig ^ are 

X3a = jj^ A; X3b = - jj^ Iq B. (36) 

Summing up all contributions and substituting them into eqn. (P), we obtain the result 
quoted in eqn. (|19|) . 
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FIGURES 

FIG. 1. The RPA series for the total transverse susceptibihty for S = infty. The first term 
represents the simple bubble, which is the building block of the ladder. 

FIG. 2. The vertex function for the exchange of magnetic fluctuations. This interaction is 
obtained from the direct transverse fermion-fermion coupling after a summation of ladder diagrams. 

FIG. 3. The 1/5 diagrams for the staggered magnetization. The solid line corresponds to a 
valence dermion, dashed line corresponds to the conduction fermion, and the wavy line represents 
the bosonic spin-wave propagator. 

FIG. 4. The diagram with two direct fermion-fermion interaction lines in both the transverse 
and longitudinal channels. 

FIG. 5. Diagrams for the irreducible susceptibility which do not contain a fermion-fermion 
interaction in the longitudinal channel. A half of the diagrams is presented. The rest are equivalent 
to those presented on Figsj5[]^ and can be obtained by replacing simultaneously all conduction 
electron propagators by valence ones and vice versa. 

FIG. 6. Diagrams for the irreducible susceptibility which contain one direct fermion-fermion 
interaction line. 

FIG. 7. The same as in Fig.^ but with two fermion-fermion interaction lines in the longitudinal 
channel. 

FIG. 8. Self-energy diagrams which contribute to the propagator of valence fermions to order 
1/5. 
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TABLES 



TABLE I. Quasiparticlo disporsion alcnig J,\,- + k,/ = tt. Euorgios aro in units of A. 



(vr/2,7r/2) 


(57r/8,37r/8) 


(37r/4,7r/4) 


(77r/8,7r/8) 


(7r,0) 


0.739 


0.751 


0.784 


0.824 


0.842 
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